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Abstract 


To safeguard against air pollution hazards to health and property, air pollution 
control programs have to be undertaken. For any air-quality-management program to be 
rational, validated simulation models of diffusion of air pollutants have to be developed. 
Such diffusion models are the yardsticks for furnishing quantitative evaluation of 
concentrations of pollutants in the ambient air. In the present work, a three-dimensional, 
time dependent, air quality grid model, has been developed by numerically solving the 
atmospheric diffusion equation. This model integrates the mass-conservation equation, 
which includes advection, turbulent diffusion and source emissions. A wind field 
inteipolation scheme is used to generate the mean transport winds from available 
measurements. While developing the scheme for numerical solution of the advection- 
diffusion equation, flux-limiter has been used so that the numerical advective flux 
function has the smoothing capabilities of the lower order scheme when it is necessary 
and the accuracy of the higher order scheme when it is possible. With the specified input 
parameters (wind speed, direction, mixing height, stability class, eddy diffusivities and 
source strength) the model is capable of predicting the gaseous pollutant concentration at 
any receptor point in the three-dimensional, space. Also, it takes into account the 
emission source located/anywhere in three-dimensional space. 

The model was validated for PMio concentrations. The model generally under 
predicted the PMio levels, primarily due to incomplete emission source inventory of 
PMio. The developed model was applied to various receptors on road-railway crossings 
in Kanpur city. Results indicate that concentration of CO and PMio can peak 2-3 times 
higher than the national Ambient Air Quality Standards at some railway crossing 
receptors due to traffic congestions. 
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1. Introduction 


1.1 Air Quality Modelling 

To safeguard against air pollution hazards to health and property, air pollution control 
programs have to be undertaken. The concentration of a pollutant at the points of 
significance from community point of view, depends on source strength and the dispersal 
capacity of the intervening atmosphere. Hence, for any air-quality-management program 
to be rational, simulation models of diffusion of air pollutants have to be developed. 
Such diffusion models are the yardsticks for furnishing quantitative evaluation of 
concentrations of pollutants in the ambient air. 

Mathematical diffusion models can be used as an aid to decision making in air 
pollution control and for city and regional planning. These models can be used in long- 
range air-source management programs for the evaluation of future air quality with 
regard to increased population. They can be used to evaluate the effectiveness of 
individual control or abatement procedures for a given source or a combination of 
sources. Models may also be useful for determining most suitable locations for pollution- 
causing industrial units as also that of sampling stations used to monitor air pollution 
levels. It may be possible to use short-range predictions of meteorological variables with 
appropriate source information to produce maps of expected air pollution concentrations 
in a few minutes with the help of a digital computer. These would be an aid in initiating 
immediate measures necessary to achieve desired air quality levels. 

A central practical problem is the prediction of the average concentration of a 
pollutant, at a specified location, which results from the contributing sources whose 
emission rates and locations are prescribed. This problem is solved only after a variety of 
simplifying assumptions are made. The understanding of the limitations that are imposed 
by these simplifications is necessary to evaluate the circumstances under which the 
various methods apply. 
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The basis of most mathematical descriptions of atmospheric diffusion is the “K- theory” 
equation, which relates the concentration of a substance, as a function of position and 
time, with the action of average wind flows and turbulence. 

The K- theory equation in the absence of chemical reaction and source terms is 


dc d(»c) 3(vc) d{wz) _ d 
dt dx dy dz dx 
Here, 



dc ' 

dx j 



yy 


dc ' 


dz 


K„ 



( 1 . 1 ) 


c 

U,v,w 

K xx , K yy and K a 

x,y,z 

t 


concentration of the pollutant 

wind speeds in the x,y and z directions respectively 

coefficients, known as eddy diffusivities 

distances along the Cartesian axes 

time 


The terms on the right hand side in Equation (1.1) represent diffusion due to 
turbulence; while on the left hand side (except the first term) represent the advection. 

To obtain the classic Gaussian plume model (analytical solution of Equation 1.1), 
the following assumptions are made: 

• The solution is time invariant 

• The wind speed is not a function of position 

• The diffusivities are not functions of position 

• Diffusion in the x-direction is insignificant compared to mean flow in that direction. 

( 

With these assumptions, the equation has an analytical solution, which expresses the 
concentration as a Gaussian or Normal distribution in the y and z directions. 

While the Gaussian equations have been widely used for atmospheric diffusion 
calculations, the lack of ability to include changes in wind speed with height and 
nonlinear chemical reactions limits the situations in which they may be used. 
Applications of simple Gaussian models, under complex urban conditions have been 
questioned by several researchers (Budianski, 1980). The atmospheric diffusion equation 
provides a more general approach to atmospheric diffusion calculations than the 
Gaussian models, since the Gaussian models have been shown to be special cases of that 
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equation when the wind speed is uniform and the eddy diffusivities are constant. To 
overcome the limitations of application of Gaussian model in complex situations, the 
above described advection- diffusion equation has been solved numerically (Numerical 
Models). 

The key problem in the use of Equation (1.1) is to choose the functional forms of 
the wind speeds, u, v, w and the eddy diffusivities, K xx , Kyy, K^, for the particular 
situation of interest. 

Appendix A shows available methods and formulations to estimate variations in 
wind speeds and eddy diffusivities with height and stability conditions. 

Here, Eulerian formulation approach is being used for the numerical solution of 
the Equation (1.1). In this formulation, a fixed coordinate system, or grid, is laid out over 
the entire region of interest. The concentration of a pollutant in each cube of the grid is 
then calculated by solving the equation, using numerical methods and with the help of a 
digital computer, over a series of small time. In the grids, which contain sources, a 
source term is added to the right side of the K- theory equation. This approach is most 
useful for situations in which there are multiple sources or for which predicted 
concentrations are needed for the entire region. 

1.2 Thesis Objectives 

In view of the above discussions concerning application of Gaussian models under 
complex urban settings, the objectives of the thesis can be worded as under 

1. To develop a three-dimensional time dependent numerical model for dispersion of air 
pollutant under complex urban settings including traffic congestion/stagnation. 

2. To apply the developed model for prediction of pollutant concentration at traffic 
intersections in the city of Kanpur characterized by complex criss-cross network of 
road and railway line with multiple level crossings. 



2. Literature Review 


2.1 Vehicular Source Air Quality Modelling 

In most of the Indian cities, traffic is the most important source of air pollution, and 
pollution from major roads is also important in suburban and rural areas. Vehicular 
dispersion models are therefore essential computational tools in modem urban planning. 

Most frequently, due to their simplicity and direct applicability for estimates on a 
local scale, various versions of the Gaussian line source model have been used for 
dispersion evaluations from a road. Such models include GM (Chock, 1978), GFLSM 
(Luhar and Patil, 1989), IITLS (Goyal and Krishna, 1999) and CAR-FMI (Harkonen et 
al., 2001). 

A simple line source (GM) model was proposed by Chock (1978) to discuss the 
downwind dispersion of pollutants near the roadway. It required only simple calculations 
and was found to be potentially more accurate than the existing Gaussian models, 
especially under adverse meteorological conditions. 

Because of the infinite line source assumption, the GM model is applicable to a 
situation where the upwind segment of the road (measured from the perpendicular line 
drawn from the receptor) is at least three times the distance between the receptor and the 
road. Hence, the GM model overpredicts the concentration by a considerable amount if 
the above constraint is not satisfied (Luhar and Patil, 1989). This constraint can be 
handled by using the HIWAY-2 model (Peterson, 1980), but its performance is not as 
good as that of GM model. 

In HIWAY-2 and CALINE-4 models, the concentrations predicted by Gaussian 
line source equation are solved by a numerical procedure. This procedure divides the 
road into a series of elements, from which incremental concentrations are then computed 
and summed up. Both models allow for a finite mixing height in the computations. 

The GFLSM model (Luhar and Patil, 1989) utilizes an analytical solution of the 
Gaussian equation for a finite line source, this solution allows for any wind direction 
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with respect to the road. The analytic solution has been originally derived from Gaussian 
formulae similar to, for instance, the HIWAY-2 and CALINE-4 models, except that the 
mixing height has been assumed to be infinite. 

Luhar and Patil (1989) presented evaluation of road dispersion models against 
experimental data sets. They presented testing of the GFLSM model, including data for 
fine paiticulate matter. The particulate matter data set was produced at a location on an 
expressway in NewYork, and also contained traffic counts and meteorological 
measurements. The model clearly overpredicted the measured fine particulate matter 
concentrations. 

Another line source model, based on Gaussian formulation, has been developed 
by IIT Delhi (IITLS model, Goyal and Krishna, 1999), in which improvements have 
been made in the calculation of source strength at major intersections in a city by 
considering the vehicular speed. 

Harkonen et al. (2001) evaluated the CAR-FMI model against measurement near 
a major road. They compared the measured concentration data with the prediction of the 
road network dispersion model. They had also analyzed the difference between model 
predictions and measured data in terms of meteorological parameters. 

Industrial Source Complex Short-Term (ISCST 3) model is a regulatory model of 
USEPA (454/B-95-003a-c). It is primarily used for point sources of industrial origin. 
However, like other Gaussian models, this can only calculate steady state concentrations 
with smallest time resolution being one hour. 

All the above-discussed models are based on the Gaussian formulation. Under a 
complex urban setting, the source location varies in three-dimensional space. Also, the 
receptor point of our interest may be anywhere in the space up to a reasonable height. 
Variation of wind speed and eddy diffusivities must be taken into account for such cases. 
With the limitations in the applicability of these models under different conditions, as 
discussed earlier, development of a numerical algorithm is needed to solve the basic 
transport equation in Eulerian grid models and such a model is developed in chapter 3. 



3. Development of Numerical Air Quality Model 


3.1 Numerical Method for Solving the Advection-Diffusion Problem 


The objective of this chapter is to devise a scheme for the solution of Equation 1.1 of the 
atmospheric diffusion, which has been reproduced here for convenience. 



to predict the time dependent concentration of air pollutants in three dimension. Here c is 
the concentration of the pollutant, u and v are the horizontal wind velocities and k xx , 1% 
and kzz are the diffusivities. There will be source term and initial and boundary 
conditions with the problem. 


3.1.1 Advection Scheme 


To analyse the various numerical advection schemes, we first consider one-dimensional 
wave equation 

§ + A(«e)- 0 (3.2) 

ot OX 

Consider grid points k-1, k and k+1. k-1/2 and k+1/2 denote the positions between them 
respectively. 


k-1/2 k+1/2 



Integrating equation (3,2) between Xk. 1 / 2 , Xk+ 1/2 and t n , t n +.i=tn + At gives 
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(3.3) 

where R=At/Ax. Here F denotes the flux through the respective faces. In the case of 
upwind scheme, F kiU2 is replaced by 


h 


’ n 

Am 1/2 




(3.4) 


This means if u>0, then /'*" 1/2 = u c" and if u<0 ,F£ U2 = u c n M . Though the scheme has 

nice property, it is a lower order scheme since it is only first order in space. Thus the 
scheme sutlers from excessive numerical diffusion. 

Hence high-resolution schemes should be applied here. High-resolution schemes are 
generally second-order accurate on smooth sections of the solution, nonoscillatory and 
capable of resolving discontinuities in the solution. Flux-Limiter methods are being 
applied for this purpose. The numerical flux function is written as 


//" 

' a »i/2 


/'■ " 4- (/>" F " 

1 I ln\l 2 ’ r*tl/2 L' Itki 


1/2 


1? n 

1 Lkk\l 2. 


(3.5) 


where Fj, and Fn are numerical flux functions of a lower order scheme and a higher 
order scheme respectively and </> ktU2 is a coefficient (called as flux-limiter) that is yet to 

be defined. It is to be noted that the flux in upwind scheme given in equation (3.4) is 
lower order. 'Hie idea is to develop a scheme where the numerical flux function has the 
smoothing capabilities of the lower order scheme when it is necessary and the accuracy 
of the higher order scheme when it is possible. Thus^ is being defined such that it is 

near 1 on smooth sections of the solution (/^" 1/2 will be approximately equal to 

h' n " kiU1 ) and near 0 on sections of the solution that contain steep gradient or 

discontinuities (/'j",, 2 will be approximately equal to F t ” ki{n ). 
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And, Equation (3.5) can also be written as 


Pk+m = Fffk+l/2 ~ (* “ ^*+1/2 )i^ r //*+l/2 ~ ^£*+1/2.1 (3-6) 

Thus the flux limiter methods can be viewed either as perturbing the flux of a lower 
order scheme or perturbing the flux of a higher order scheme. 

We already know that upwind flux can be taken as F L . For the higher order 
flux, following schemes can be used successfully. 

The first one is called Lax- Wendroff scheme: 

' r „;.„ 2 =“c;+L‘( (3.7) 

where 8,c n k = c n M -c n k . 

The second scheme is called Mid-Point scheme where the flux is gi/en by 

+»;.,) (3-8) 

There are many ways to accomplish the goal of achieving <j> near 1 in smooth 
regions and 0 near discontinuities. A common approach is to write 

=<*(»;); 

where 0 k is a smoothness parameter. We take 


n 




(3-9) 


where d\c k = c k - c k „ , . 

6 k is used to measure the smoothness of the solution, because when the solution is 
not changing much, it will be approximately 1, and when the solution is changing a 
significant amount, it will be very different from l.The function <f> is chosen such 
that the numerical scheme satisfy some desirable properties. These will be satisfied 

if we require $ to satisfy 


0 :£ 


m 


£2 


and 0<,</>{d)<2 


& 


(3.10) 
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for all 9. Some of the choices for limiter functions that are in the literature include 
the following: 

Superbee Limiter: 

^(#) = max {0, min [l, 2d\ min [<9, 2]} (3.1 1) 


Van Leer Limiter: 




0 + | 0 | 
1 + 1*1 


(3.12) 


The above numerical scheme can easily be extended to two dimensions. In this case, 
wc consider 


~c + ~—{uc) + -~(vc) = 0 (3.13) 

dt dx ’ dy x ! v ' 

Integrating the above equation with respect to t from t n to t n +i, from xj_j a to Xj+ 1/2 

with respect to x and from yv..\n to yk+1/2 with respect to y, we get 


Cm = <• -R'[ [Q- M „ n -B m-„A ( 3 . 14 ) 

Where R x = A(/Ax and R v - At/Ay. Here P and Q are the fluxes in the x and y 
directions respectively. Proceeding as before, we write 


Pj*U2,k ~ )*M2.k + fijun.k 


p n __ p n 

r H J+lf2,k r Lj+\n,k. 


where <f> is the limiter function. Pi, and Ph denote the fluxes based on lower and 
higher order schemes. Similarly we write the flux along the y- direction also. The 
same scheme is obviously generalized to three dimensions. 



3.1.2 Diffusion Scheme 


Again, first we consider a simple one-dimensional equation of the form 


— - — f K 

dt dx v dxj 


(3.15) 


Where K is the diffusivity constant. Integrating the above equation from t n to t„+i /2 
for t and from Xj. 1/2 to Xj+ 1/2 for x, we get 


Ax 

~At 




dc 

dx 


J+U2 


-K 


J- 1/2 


dc 

dx 


j- 1/2 


Expanding the term on the right hand side further, we get 


(3.16) 


Ax r 
At 


i c j 




$]■ 


K 


*7+1 






;+i/ 2 


Ax 


K 


„»+■ 

V-1 


J-U2 


Ax 


(3.17) 


This is an implicit scheme, since the time step for the explicit scheme is very small. 
The above equation can be written as 


a J c T 


+«>♦> +a H c *:i 


where 


a /+U2 

a j-il2 


~ Kj+U2 
= Kj-U2 


_At_ 

Ax 2 

A? 

Ax 2 


and 


(3.18) 


a j - 1 + a }+\n + a j-\n 

The above technique can be generalized to two and three dimensions in obvious 
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3.1.3 Solution Method 


The Equation (3.1) in three dimensions is discretized as (using the techniques of the 
previous sections) 


®ij.k ^i.j,k ~~ ^t,j\k 


r,\gi 


>- I 
«+l 


ij+l/ 2 *k 




A 


a i+\n.j,k C i+\J % k + a i-l/2J,k C t-\J,k a iJ+U2 y k C iJ+\,k + a tj-\/2,k C ij-Uk 


+ a iJ,k+U2 


C 


n+\ 

t,j\k + 1 


+ a iJ.k~ 1/2 


C 


n hi 


(3.19) 


Where R’s and a’s have expressions similar to those of the previous sections. In the 
scheme, advection is treated explicitly and diffusion is treated implicitly. The above 
equations on the mesh are solved using iterative techniques such as Gauss -Seidel 
method with relaxation. 


3.2 Validation of the Code for Advection Algorithm 


3.2.1 Standard Cone Test 

To verify the advection scheme using a limiter, the classic cone test in a solid body 
rotational flow-field has been used. The flow-field velocities are defined over a 
1 00x1 00-grid cell domain: 

u = -©y and v = cox with © = 0.1. 

The initial condition has a cone of radius 15 grids and height of 5, above a 
background of base value 1. The cone is centered at (50,75). This test is often used 
to investigate the accuracy of advection schemes (Clappier, A., 1998). It is 
important to note that this rotational flow-field is monodimensionally divergence 

free 
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( *1 = ^ = 0 ) 
v dx dy 

During this test, the cone should rotate around the domain, keeping the same shape. 

3.2.2 Suitability of the Limiter and High Order Schemes 

To verify the best suitability of the limiter and high order schemes, the combinations of 
the earlier mentioned limiters with high order schemes has been taken one by one and 
standard cone test is applied. The combination, which preserves the shape of the cone at 
most with less reduction in peak is taken to be the best and is applied here. 

Here, we find that the Van-leer limiter with lax-Wendioff scheme fulfills the 
requisites better in comparison to the others. The cone test results for one full rotation in 
the horizontal plane, applied to the various combinations are as given here: 


1-1; one rotation 


"ap81 jdat" 



Figure 3.1 Mid-point Scheme; Van-leer Limiter 



'ap81 .dat* 



Figure 3,2 Lax-Wcndroff Scheme; Van-leer Limiter 




"apsixlat" — 



Figure 3 3 Mid-point Scheme; Superbee Limiter 
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2,2 


"ap81.dat" 



Figure 3.4 Lax-Wendroff Scheme; Superbee Limiter 


3.2.3 Cone Test for the Lower Order Scheme 

As discussed earlier, lower order scheme is the Upwind Scheme. When cone test is 
applied for the scheme, there is much reduction in the peak with deformation in the 
basic shape of the cone. 

The shape is no more preserved after one and two rotations of the cone. 

Also, the peak of the cone has reduced considerably (sixty percent reduction after 
one rotation and seventy five percent reduction after two rotations), thus restricting 
for the use of the lower order scheme. 
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Upwind Sc hams; No Limiter, two rotations 



Figure 3.7 After two rotations of the cone in horizontal plane 


3.2.4 Cone Test for the Limiter used in the Model 

While using the combination of Van-leer limiter and Lax-Wendroff scheme, we find 
the following results after single and two rotations of the cone in the horizontal 
plane. As in the previous case, the cone has been introduced with a peak of 5 units 
over a base value of 1 unit. After rotation in the horizontal plane, the peak goes to 
4.5 (i.e. 90% of the original) showing that little numerical diffusion is taking place. 
And, after two full rotations in the horizontal plane, it goes to 4.0 (i.e. 80% of the 
original). Furthermore, after two rotations of the cone, the peak of the cone is not 
going to decrease significantly, indicating that afterwards, the decrease in the peak 
concentration is less, and thus making the scheme to be acceptable. 




Figure 3.8 initial cone (Cen 
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lax-W»nd lofl Sch«mo, Van-lw Umiter;two rotations 



Figure 3.10 After two rotations in the horizontal plane 


4. Parameters Assignment and Application of Model 


This chapter presents the description of various parameters used in model application 
and details of field sampling . 

4.1 Meteorological Parameters 

4.1.1 Wind Speed and Direction 

The wind speed and direction data were obtained from Weather Monitoring Station of 
Chandra Shekhar Azad Agricultural University (CSAAU), Kanpur and Weather 
Monitoring Station maintained at Environmental Engineering (EE) Laboratory, IIT 
Kanpur. 

Measurement of wind speeds is done using a three-cup anemometer. The 
anemometer has a vertical axis and three cups, which captures the wind. The number of 
revolutions per minute is registered electronically. The anemometer is fitted with a wind 
vane to detect the wind direction (Envirotech Model WM-251). 

The wind data are in the form of 24 hr average value at CSAAU and online 
continuous data of wind are available from EE laboratory, IIT Kanpur. Wind data for 
specific days, for which field sampling was done, were collected from these two sources. 

4.1.2 Atmospheric Stability Class 

Atmospheric stability classes were obtained using the wind speeds and climatic 
conditions (such as cloud cover, solar radiation etc.) using the insolation based technique 
(Seinfield, 1998). Climatic condition data were obtained from CSAAU Weather 
Monitoring Station. 
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4.1.3 Eddy Diffusivities 

Eddy diffusivity in vertical direction is calculated using the wind speed and mixing 
height values (for unstable conditions) as per the methods presented in Appendix A. 

In horizontal direction, diffusivities are taken as constant, equal in both (along the 
x & y) directions. Its value is taken equal to 10 m 2 /s (Seinfeld, 1998). 

4.2 Source Strength 

4.2.1 Emission Sources 
Vehicles 

The primary interest of application of the model is in prediction of pollutant 
concentrations mainly from vehicular emissions. Emission rates have been estimated 
using number of vehicles and emission factors developed by Indian Institute of 
Petroleum (IIP). (Table 4.1.; IIP, 1982; taken from Sharma, 1999). 

Sivacoumar and Thanasekaran, (2000) have proposed emission factors as a 
function of speed of vehicle and the same are used in calculating the emission quantities 
of carbon monoxide, when vehicles are in motion. 

Traffic pattern for the different types of vehicles have been shown in Figures 4. 1 
and 4.2 at two locations, Gurdev crossing and Rawatpur crossing (Figure 4.3).It is clear 
from these figures that traffic density is much higher at Rawatpur than at Gurdev 
crossing. There are generally two peaks, in morning between 10-11 AM and in evening 
between 6-7 PM. Vehicle pattern shows that two wheelers are maximum in number 
followed by Vikram Tempos and four wheelers.This information on number of vehicles 
is translated into hourly emission rates at these two crossings. 

It may be noted that at each of the above crossings, two receptors points have 
been taken to predict air concentrations (GD-R1 and GD-R2 at Gurdev crossing and RP- 
R! and RP-R2 at Rawatpur crossing; Figures 4.4 and 4.5). 



Table 4.1 Emission factors for different vehicles (IIP, 1982) 


Vehicle Type 

CO 

HC 

NO* 

TSP* 

Two- wheelers 

8.30 

5.18 

— 

— 

Cars 

24.03 

3.57 

1.57 

— 

Three-wheelers 

12.25 

7.65 

— 

— 

Buses (Suburban) 

5.87 

2.27 

11.10 

0.37 

Buses (Urban) 

4.51 

1.75 

8.52 

0.28 

Trucks 

3.52 

1.36 

6.66 

0.22 

Light Commercial 

1.30 

0.50 

2.50 

0.10 


Vehicles 


* : TSP has been taken as PMio (Particulate matter less than 10 pm size) 
N.B. All the values are in g/km/vehicle. 
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Figure 4.1 Traffic Pattern at Gurdev Crossing 
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Traffic pattern at Rawatpur 
Trucks «* Buses 
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Traffic Pattern at Rawatpur: Four Wheelers 
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Figure 4.2 Traffic Pattern at Rawatpur Crossing 









Figure 4.3 Map of Kanpur City showing the measurement site locations 
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Figure 4,4 Road Layout at Gurdev Crossing with receptor locations 
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To Gurdev 



From Gol 
Chouraha 



Figure 4.5 Road Layout at Rawatpur Crossing with receptor locations 
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4.2.2 Road Side Dust Emission 

One major source of fine particulate is road dust resulting from movement of traffic. 
USEPA has developed generic equations to estimate particulate emissions from paved 
and unpaved roads depending on road conditions and weight of vehicle. 

PM io and PM2.5 Emission Factor Equations 

Dust emissions from paved roads have been found to vary with the silt loading present 
on the road surface as well as the weight of vehicles traveling on the road. For unpaved 
road, besides the silt content and vehicle weight, the dust emissions also vary with 
vehicle speed and average number of vehicle wheels. 

Paved Road 

Road dust emission due to traffic has been calculated using the following empirical 
equation (EPA-website) 

E = k* (SL/2) 0 ' 65 (W/3) IS ((365-N)/365) 
where E = particulate emission, g/vkmt 
k = 4.6 for PM 10 , 1.1 forPMzs 
SL * silt loading, g/m 2 (measured in field as 15 gm/m 2 ) 

W = average vehicle weight, tons (taken as per actual weights) 

N = number of days in year with having rain more than 2.5 cm. (taken as 40) 

The derived emission factors for road dust emission corresponding to paved 


roads are shown in Table 4.2. 



28 


Table 4.2 Fugitive Emission Factors for road-side dust for different vehicles 


Type of Vehicles 

Fugitive Emission Factor (g/km) 

Truck 

34.825 

Bus 

26.492 

Vikram 

4.353 

Two-wheeler 

0.074 

Three-wheeler 

0.599 

Four-wheeler 

6.702 


4.2.3 Background Pollutant Concentration 

Background concentration of PMio in atmosphere was taken as 60 pg/m 3 as estimated by 
S harm a et al. (2002) for Kanpur city. 

4.3 Measurement Sites and Field Sampling 

4.3.1 Site locations 

The emission inventory of PMio in the city of Kanpur has not been prepared. The main 
sources of pollutant emissions in the study area include industrial, vehicular and other 
combustion processes. The air pollutant emissions are from a coal-based 200 MW 
thermal power plant (industries) and from 250000 vehicles, which ply on the road. 

Figure 4.4 shows receptor locations (GD-R1 and GD-R2), where actual sampling 
of PM io was undertaken during the month of January 2002 using high volume sampler 
(Ncttlc,Mumbai*NPM~HVS-850). 





4.4 Pollutant Parameters 


In this study, only two pollutant parameters are considered for prediction of ambient air 
concentrations. These are: 

(i) PMio;and 

(ii) Carbon monoxide (CO) 

Reasons for selecting these to pollutants include that PMio levels are very high in 
the city of Kanpur (about 2-4 times higher than National Ambient Air Quality Standards) 

and CO is the main pollutant from automobiles. 

4.5 Grid Description for Model Application 

To solve equation 3.1 numerically, a three dimensional space with 
dimensions 500 m x 500 m x 200 m was taken and was divided in the grids 50 x 50 x 40. 
Thus the size of the small grid is 10 m x 10 m x 5 m. 



5. Results and Discussion 


5.1. Model Validation 

As stated in Chapter 4 (Section 4.3), two locations (GD-R1 and GD-R2) were taken for 
the sampling of PM to to generate the data for model validation. Table 5.1 presents the 
results of measured PM to levels at these two locations. 


Table 5.1 : Measured PM to levels (8- hr average) at GD- R1 and GD- R2 


Date 

PMio (fig/m 3 ), GD- R1 
(8-hr average) 

PMio (pg/m 3 ), GD- R2 
(8-hr average) 

January 25, 2002 

428 

98 

January 28, 2002 

477 

149 

January 30, 2002 

617 

118 

January 31,2002 

354 

145 


It is clear from Table 5.1 that as one moves away from the road, the levels of 
PM to drop rapidly. It may be observed that measured levels of PMio at both locations 
exceed the National Ambient Air Quality Standard of 100 (ig/m 3 . 

ft may be mentioned that it is not desirable to take a sample smaller than 8-hr 
duration, as accuracy of the gravimetric analysis reduces for shorter sampling. 
Therefore, model validation can only be done for 8-hr measured values vis-^-vis 8-hr 
average model-computed values. It may, however, be noted that model has the 
capability of predicting concentrations for very small time steps (in seconds), if desired. 

Figure 5.1 shows the linear plot between measured and model-computed 
8-hourly values at the two receptors considered for model validation. 
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FMto Concentration (pg/cu*m); model-computed 


Figure 5.1 Comparison of Measured and Predicted PMio Concentrations 


There ate several m'loworih) points that can be drawn from Figure 5.1. 

• The model consistently under predicts the concentration compared to the measured 
comcmi atoms In fact, this pattern of under prediction was expected. It may be noted 
that the model uses the emissions caused by vehicles and road dust alone, whereas 
the measured loncenlrafion is the resultant of the combination of emissions from 
several soutces in the city of Kanpur. 

• Figure 5 I suggests that even if predicted value is zero (i.e. no movement of 
vehicles), a kuTromnd level of 71.33 jig/m 1 will exist at these two locations, which 
can be attributed to contribution from the other sources. 

• rite model prediction and actual pattern of PMio concentration match very well as 
indicated by R* value ( T 0 %) and slope of the line (=1.0972) being close to unity, 
The value of R* is significant at 5% level of significance in a statistical sense. 
Therefore, it can he concluded that model adequately describes the dispersion of 
putlul.mk in urban environment 
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5.2 Model Application 

5.2.1 Results at GD- R1 and GD-R2 

Figures 5.2 and 5.3 present hourly variation in the value of PMio at two receptors (GD- 
R1 and GD- R2) at Gurdev crossing. As expected, the predicted concentrations are 
higher at location GD- R1 than that at GD- R2 on both the days of prediction and this 
situation is true for every hour. It suggests that the impact of vehicular emission 
diminishes, as one moves av'ay from the roadside. This pattern of decline in air 
concentration has also been found by other researchers (Harrison and Johnston, 1985; 
Hewitt and Rashed, 1991 ) for PAHs and Pb (emitted from vehicles). 

It can be noted that peak hourly concentrations can be several times higher than 
National Ambient Air Quality Standards (i.e. 100 pg/m 3 ), which may cause discomfort 
and other health impacts to the people those are exposed for a long time to the air on the 
roadside. 

The model results of PMio concentrations follow the vehicular density for a given 
hour; more vehicles at the intersection, higher is the concentration. 

Another interesting fact, which emerges out from the comparison of Figures 5.2 
and 5.3 is that PMio level shooted up on January 28, 2002 compared to January 25, 2002 
(in spite of same emission rate). To understand this variation, one need to examine the 
meteorological parameters such as values of average wind speed (U) and A0 (difference 
in wind direction and the direction of receptor, in degrees). As the wind speed increases, 
the concentration comes down due to dilution effect and as A0 becomes small, it suggests 
that receptor is coming in the downwind direction of emissions and will receive more 
pollution. In case of January 28, 2002 results for GD- Rl, U and A0, both have decreased 
(compared to values on January 25, 2002) resulting in high level ofPMjo at GD- Rl. 

The influence of meteorology can also be seen at GD- R2 to a lesser extent, while 
comparing the concentrations on January 25 and 28, 2002. 



PM10 Concentration fmicrogram/cu-m) PM10 Concentration (microgram/cu-m) 
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Figure 5.2 Diurnal variation in PM t0 Concentration as predicted by the Numerical 
Model (25th Jan 2002) 

(tl Average wind speed, AO -> Dir, of wind - dir. of receptor) 
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Figure 5.3 Diurnal variation in PM l0 Concentration as predicted by the Numerical 
Model. (28th Jan 2002) 

(U -> Average wind speed, AO -> Dir. of wind - dir. of receptor) 
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5 2.2 Results at RP-R1 and RP-R2 

At these two locations of Rawatpur crossing, the traffic density is much higher than at 
Gurdev crossing (see Figures 4.1 and 4.2). At this crossing also, two receptor locations 
RP- RI (close to the intersection) and RP- R2 (away from the intersection) are taken for 
prediction of concentration. Model predictions of concentrations of two pollutants, CO 
and PMio have been carried out for April 10, 2002. 

Attempts have also been made to compare the results of the developed model 
with that of ISCST 3 model of USEPA (USPEA,1995). Figures 5.4 and 5.5 show the 
comparisons made between the predicted values of Carbon monoxide at receptor 
locations RP- Rl and RP- R2 using the developed numerical model and ISCST 3 model. 




Figure 5.4 Comparison of CO concentration at RP- R1 using the 
Numerical Model and ISCT 3 Model (April 10, 2002) 
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It can be observed ftom the Figures 5.4 and 5.5 that the CO concentration values 
obtained using the numerical model is comparatively higher than that of ISCST 3 model. 
The concentration values are in the range of 1-1 1 mg/m 3 .The range is comparable to 8 hr 
CO concentration in Delhi, i e 5 mg/m 3 (CPCB, 2000). The sharp peaks in the results of 
numerical model is because of the fact that the gates of level crossing were closed due to 
arr j va l 0 f trains, When gates are closed there is considerable increase in traffic (see 
Figure 4.2), which causes the consequent increase in the pollutant level. This pattern of 
sharp increase cannot be obtained using the Gaussian models (e.g. ISCST 3 model). 
Further one can get an idea that the numerical model-predicted values are much higher 
than concentrations predicted by ISCST 3 model. Since, it was found that numerical 
model mulct predicts { see section on model validation), it is reasonable to conclude that 
the developed model perform better than the ISCST 3 in the present situation with the 
advantage of estimation of short-term peak concentrations. 

5.2.3 Prediction of PM 10 Concentration at RP- Rl and RP- R2 

The numerical model pmlii ted values at the locations RP- Rl and RP- R2 fot PMio 
concentrations have been shown in Figure 5.6. 



Figure 5.6 Numerical Model Predicted PMio values at RP- Rl and RP- R2 
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It can be observed from Figure 5.6, that the PM !0 levels at the location RP- R1 is 
higher than that at RP- R2 because RP-R1 is nearer to the intersection and receives more 
pollution. 

The PM hi concentration at both the receptors exceed the National Ambient Air 
Quality Standard of 100 pg/m\ While the average concentration of PMio is 450-500 
pg/m ? , the short- term peaks may even exceed the value of 800 pg/m 3 ofPMio, which is 
matter of concern for public health. 
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6. Conclusions 


The following conclusions can be drawn from present study 

The developed numerical model can capture the short-term peaks of air pollutant 
concentration arising due to unsteady urban conditions prevalent in city of Kanpur due to 
criss-cross network of road and railways. The model results clearly show that even 
though number of vehicles in the city of Kanpur is much less compared to large city like 
Delhi, pollutant concentration in Kanpur can be much higher due to traffic congestion 
and pool planning. 

In terms of numeric concentration, model suggests that carbon monoxide 
concentration exceed the national ambient air quality standard by a factor of 2 at 
Rawatpur crossing Similarly, PMio levels also exceed the standard by a factor of 4-6 

on the two crossings studied in Kanpur city. 



7. Scope for Future Work 


There is an urgent need for good quality experimental data, which should be specifically 
produced for further model validation and evaluation purposes. 

The present mode! docs not incorporate reactions of chemical species. A more 
detailed output can be obtained after consideration is made for the formation and 
removal of pollutant species of interest in the atmosphere. 

The model should be used to develop pollution control strategies for various 
crossings in the city of Kanpur. 
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APPENDIX A 


Wind Speed Variation and Eddy Diffusivity 


Variation of Wind with Height in the Atmosphere 

The atmosphere near the surface of the earth can be divided into three layers — 
the free atmosphere, the Ekman layer, and the surface layer. The Ekman layer extends to 
a height of from 300 to 500 m depending upon the type of terrain, with the greater 
thickness corresponding to the more disturbed terrain. 

The layer immediately adjacent to the surface, typically upto 30 to 50 m from the 
ground, is called the surface layer. Within this layer, the vertical turbulent fluxes of 
momentum and heat are assumed constant with respect to height, and indeed they define 

the extent of this region. 

Considering two-dimensional turbulent flow of air in the surface layer parallel to the 
ground at a height 23, the wind speed u (Z3) can be given by, 


«(*>) = l ln £j 

U. K Z 0 


Z 3m Zq 


(A.l) 


where k is a constant proportionality factor with the experimental value of 0.4 and zo is 
called the roughness length and values of zo for typical surfaces are given in table 1.1. 


Table A.1 Roughness Lengths for Various Surfaces (Seinfeld, 1998) 


Surface 

Z 0 (m) 

Very smooth (ice, mud flats) 

10' 5 

Snow 

10' 3 

Smooth sea 

10' 3 

Level desert 

10' 3 

Lawn 

10' 2 

Uncut grass 

0.05 

Fully grown root crops 

0.1 

Tree covered 

1 

Low density residential 

2 



Of the three parameters K, z<) and u* mentioned above, zo is the property of the roughness. 
To calculate the value of u* at different depths, wind speed is measured at some 
convenient height, say z. Then equation (A.l) can be used to calculate the value of u* 
and subsequent specification of wind speeds at all values of heights. 

The equation (A.l) corresponds to neutral weather conditions. For stable and 
unstable weather conditions, we can use the equation 

«(z 3 ) _ 
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In 


f 7 


Inj 


V *o J 
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4.7 i , 
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stable 


\ z o J 


+ In 


k + l )k + O 2 

(Vr+ l )k + l ) 2 


+ 2 (tan ?] r - tan 7 ] a ) unstable 


(A.2) 


where r) r and t]o are given by 


I 

tA 

l 

1/4 

and 

% = 

1 - 15 — 

L 



L 


h, called Monin-Obukhov length is not a parameter that is routinely measured. It is the 
height above the ground at which the production of turbulence by both mechanical and 
buoyancy forces is equal. It provides a measure of the stability of the surface layer. To 
fulfil the need, we can approximate the ({older’ s (1972) plot between the Pasquill 
stability classes, the roughness length zq and L, by the correlation 

L = (a + b log zo)' 1 (A.3) 

The values of a and b have been listed in Table A.2. 

Determination of Pasquill stability classes is made as per Table A.3. 
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Table A.2 Coefficients for straight-line approximation to Golder’s plots as a 
function of stability classes (Seinfeld, 1998) 



Pasquill Stability Class 

Coefficients 

a b 

Extremely unstable 


A 

-0.096 

0.029 

Moderately unstable 


B 

-0.037 

0.029 

Slightly unstable 


C 

-0.002 

0.018 

Neutral 


D 

0 

0 

Slightly stable 


E 

+0.004 

-0.018 

Moderately stable 


F 

+0.035 

-0.036 

Table A.3 Estimation of Pasquill Stability Classes (Seinfeld, 1998) 


Surface Wind 

Speed at 10m 

Solar Radiation 

Night time 

fraction 

cloud Cover 

(m/s) Strong 

Moderate 

Slight 

" >= 4/8 

<=3/8 

<2 A 

A-B 

B 

— 

— 

2-3 A-B 

B 

C 

E 

F 

3-5 B 

B-C 

C 

D 

E 

5-6 C 

C-D 

D 

D 

D 

>6 C 

D 

D 

D 

D 
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Eddy Diffusivity Coefficients 

Horizontal Eddy Diffusivity Coefficient K xx and K w : 

Equation (1.1.1) shows that the turbulent diffusivity has three components, but 
the more caielul formulation demonstrates that it is a second order tensor with nine 
components. Thus, for example, the x- component of the Reynolds dispersion terms are 

more rigorously given by 


A, 


f , A' — - + K ~ 

" rip " f)z 

However, no information on the otT-diagonal terms, such as K xy and K xz are available, 
and it is customary to assume that they are zero. Furthermore, we usually assume that the 
horizontal eddy diffusivities are equal, K xx = K yy . The equation for turbulent dispersion 
then contains only K xx , K yy and K zz terms. 

For a range of typical meteorological conditions, we find diffusivities under 
unstable conditions in the order of (50-100 m 2 /s). 


Vertical Eddy diffusion Coefficient K^: 

The expulsions available for K„, are based on Monin-Obukhov similarity theory 
coupled with observational or computationally generated data. It is best to organize the 
expressions according to the type of stability. 

Since we generally need expressions for K zz that extend vertically beyond the surface 
layer, we should consider some available correlations for the entire Ekman Layer. 


Unstable conditions: 

When turbulence in the atmospheric boundary layer is maintained largely by 
buoyant production, the boundary layer is said to be in a convective state. The source of 
buoyancy is the upward heat flux originating from the ground heated by solar radiation. 
Convective turbulence is relatively vigorous and causes rapid vertical mixing in the 
atmospheric boundary layer. 

In the unstable weather conditions, there is usually an inversion base height at 
z~ Zi that defines the extent of the mixed layer. The two parameters that are key in 
determining Kzz are the convective velocity scale w* and z*. An empirical 



expression for Kzz under unstable conditions, using the numerical turbulence 
model of Deardorff (1970) is as follows: 


0<— <0.05 

0.05 <—<0.6 

0 . 6 <— < 1.1 
z i 

— > 1.1 

(A.4) 

Neutral Conditions: 

The following relationship (Shir, 1973) for Kzz can be used successfully under 

neutral conditions: 


K u = k k, z exp 

where, f is the coriolis parameter. 

Stable Conditions: 

Under stable conditions, the appropriate characteristic vertical length-scale is L. 
Businger and Arya (1974) proposed a modification of surface layer similarity 

theory to extend its vertical range of applicability. 

,, ku.z ( Bfz 

& 0.74+ 4.7(2/ L) u. 






APPENDIX B 


Table Bl: Traffic Count at Gurdev Crossing (along G.T. road) 


Time 

Trucks 

Buses 

Vikrams 

2-Wheelers 

3-Wheelers 

4-Wheelers 

Tractors 

















8 a.m. - 9 a.m. 

80 

26 

352 

682 

9 

161 

3 

9 a.m. - 10 a.m. 

48 

14 

315 

778 

23 

217 

8 

10 a.m. - 11 a.m. 

46 

20 

323 

802 

21 

270 

9 

11 d.m - 12 noon 

145 

39 

321 

737 

17 

224 

12 

. 12 _nocnj.!.p- m - 

140 . 

52 

362 

712 

15 

267 

11 

1 p.m. - 2p.m. 

105 ’ 

38 

333 

678 

21 

223 

12 

2 p.m. - 3 p.m. 

113 " 

26 

289 

582 

36 

212 

10 


118 

34 

256 

537 

32 

218 

6 

4 p.m. - 5 p.m. 

95 

23 

346 

803 

28 

223 

8 

5 p.m. - 6 p.m. 

57 

15 

331 

793 

41 

231 

10 

6 p.m. - 7 p.m. 

49 

14 _ ___ 

289 

741 

35 

188 

7 

7 p.m. - 8 p.m. 

§2 , 

19~ 

322 

677 

16 

338 

12 


Table B2: Traffic Count at Rawatpur Crossing (along G.T. road) 


Time 

Trucks 

Buses 

Vikrams 

2-Wheelers 

3-Wheelers 

4-Wheelers 

Tractors 





■ ■ ' 1 

. . 


(Cars etc.) 










8 a.m. - 9 a.m, 

74 

21 

398 

910 

13 

234 

4 

9 a.m. - 10 a.m. 

42 

. • 12 - 

342 

1006 

38 

318 

8 

10 a.m. - 11 a.m. 

35 

17 

359 

1096 

37 

347 

7 

11 a.m. - 12 noon 

113 

43 

353 

1117 

28 

312 

15 

12 noon - 1 p.m. 
Ijp.m. - 2p.m. i 
2 p.m. -3 p.m. 1 

115 

48 V„ 

402; 

997 

32 

335 

8 

| ' , 88 ’ J 

34“ 

367 

934 

39 

316 

15 

I 96 

32 

325 

915 

43 

298 


3 p.m, - 4 p.rn. 

i 86 

46 

290 

884 

52 

295 

| 

4 p.m. - 5 p.m. | 

i 83 


436 

1084 

46 

321 

HUH 

5 p.m. - 6 p.m. | 

6 p.m. - 7 n.m. 1 

1 . - ■ ~i 

j 45 ] 


387 

1113 

61 

292 

8 

”40 

13 

386 

1046 

31 

335 

7 

7 p.m.- 8 p.m. | 48 

“ 14 

378 

942 

18 

370 

11 
































































Table B3: Traffic Count at Rawatpur Crossing (on the road orthogonal to G.T, road and towards 
Kakadev) 
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Train Arrival Time 




£ 

CO 

to 

CM 

00 

9:35 a.m. 

10:05 a.m. 

11:10 a.m. 

12:35 D.m. 

£ 

Q 

O 

T— 

2:45 D.m. 

E 

Q 

O 

CM 

CO 

4:40 D.m. 

5:15 p.m. and 5:45 D.m. 

6:33 p.m. 

E 

d 

o 

! 4- Wheelers 1 Tractors 




CM 

r~ 

CO 

to 

CO 


CO 

CM 

h- 

22 


CD 


(Cars etc.) 



245 (21) 

267 (58) 

214(34) 

314(27) 

308 (35) 

279 (43) 

268 (32) 

231 (29) 

257 (45) 

298 (34)(28) 

CD 

CO 

M- 

O 

CO 

288 (32) 

1 

1 

* 

e 

% 

00 




CO 

CM 

45 (16) 

65(18) 
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48(6) 

22(4) 

s 

co 

CD 

32 (1) 

CM 

CO 

to 

2 

00 

T— 




CM 

s 

CM 

to 

518 (38) 

556 (49) 

456 (32) 

423 (19) 

487(17) 

378 (24) 

627(55) 

I 

CO 

CO 

o 

s 

CO 

s 

467(36) 

g 

<0 

156(13) 

114(10) 

s 

s 

84(23) 

34(21) 

23 (13) 

to 

T** 

8 

112(30) 

56 (6)(16) 

106(12) 
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CQ 

<u 

1 



- 


V* 

S3 

CO 

t** 

N 

CO 

CM 

r* 

CO 

r~ 

CO 

CO 

s 

CO 

o 

CM 

CO 

to 


M- 

£ 

rd 

05 

k 

d 

CO 

o 

£ 

cd 

o 

T— 

£ 

d 

o> 

to 

CM 

N 

to 

*«fr 

N 

CO 

00 

CO 


10a.m. - 11 a.m.' 

C 

o 

o 

c 

CM 

T- 

E 

<0 

T m 

12 noon - 1 p.m. 

£ 

d 

CM 

t 

£ 

d 

r* 

2 p.m. - 3 p.m. 

3 p.m. - 4 p.m. 

4 p.m. - 5 p.m. 

£ 

d 

CO 

1 

£ 

d 

LO 

6 p.m. - 7 p.m. 

£ 

d 

00 

k 

d 


N.B. Values in the parentheses show the vehicles waiting on either side of the level crossings at the time of arrival of the train. 




Table B4: Traffic Count at Gurdev Crossing (on the road orthogonal to G.T, road and towards Zoo) 
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